*
* $Id$
*
* $Log: gnslwd.F,v $
* Revision 1.1.1.1  2002/06/16 15:18:40  hristov
* Separate distribution  of Geant3
*
* Revision 1.1.1.1  1999/05/18 15:55:20  fca
* AliRoot sources
*
* Revision 1.1.1.1  1995/10/24 10:21:15  cernlib
* Geant
*
*
#include "geant321/pilot.h"
*CMZ :  3.21/02 29/03/94  15.41.38  by  S.Giani
*-- Author :
      SUBROUTINE GNSLWD(NUCFLG,INT,NFL,TEKLOW)
C
C *** NEUTRON TRACKING ROUTINE FOR ENERGIES BELOW THE CUT-OFF. ***
C *** TAKE ONLY ELASTIC SCATTERING, NEUTRON CAPTURE            ***
C *** AND NUCLEAR FISSION.                                     ***
C *** NVE 11-MAY-1988 CERN GENEVA ***
C
C CALLED BY : GHEISH
C ORIGIN : H.FESEFELDT (ROUTINE NSLDOW 20-OCT-1987)
C
#include "geant321/gctrak.inc"
C --- GHEISHA COMMONS ---
#include "geant321/mxgkgh.inc"
#include "geant321/s_consts.inc"
#include "geant321/s_curpar.inc"
#include "geant321/s_result.inc"
#include "geant321/s_blankp.inc"
#include "geant321/s_prntfl.inc"
      DIMENSION RNDM(2)
C
C --- FLAGS TO INDICATE THE NUCREC ACTION ---
C NUCFLG = 0 ==> NO ACTION BY NUCREC
C          1 ==> ACTION BY NUCREC ==> SPECIAL TREATMENT IN GHEISH
      NOPT=0
      NUCFLG=0
C
C --- IN ORDER TO AVOID TROUBLES CAUSED BY ARITHMETIC INCERTAINTIES, ---
C --- RECALCULATE SOME QUANTITIES. TAKE KINETIC ENERGY EK AS MOST ---
C --- RELEVANT QUANTITY. ---
C
C --- VERY LOW KINETIC ENERGY ==> NEUTRON CAPTURE ---
      IF (EK .LT. 1.E-9) GO TO 22
C
      EN=EK+ABS(AMAS)
      P=SQRT(ABS(EN*EN-AMAS*AMAS))
      PU=SQRT(PX**2+PY**2+PZ**2)
      IF (PU .GE. 1.E-9) GO TO 7
C
      PX=0.0
      PY=0.0
      PZ=0.0
      GO TO 22
C
 7    CONTINUE
      PX=PX/PU
      PY=PY/PU
      PZ=PZ/PU
C
C --- SELECT PROCESS ACCORDING TO "INT" ---
      GO TO (23,23,21,22), INT
C
C *** NUCLEAR FISSION ***
 21   CONTINUE
      ISTOP=1
      TKIN=FISSIO(EK)
      GO TO 9999
C
C *** NEUTRON CAPTURE ***
 22   CONTINUE
      ISTOP=1
      CALL CAPTUR(NOPT)
      GO TO 9999
C
C *** ELASTIC AND INELASTIC SCATTERING ***
 23   CONTINUE
      PV( 1,MXGKPV)=P*PX
      PV( 2,MXGKPV)=P*PY
      PV( 3,MXGKPV)=P*PZ
      PV( 4,MXGKPV)=EN
      PV( 5,MXGKPV)=AMAS
      PV( 6,MXGKPV)=NCH
      PV( 7,MXGKPV)=TOF
      PV( 8,MXGKPV)=IPART
      PV( 9,MXGKPV)=0.0
      PV(10,MXGKPV)=USERW
C
C --- SPECIAL TREATMENT FOR INELASTIC SCATTERING IN HEAVY MEDIA ---
      IF ((INT .EQ. 2) .AND. (ATNO2 .GE. 1.5)) GO TO 29
C
C *** ELASTIC SCATTERING ***
 30   CONTINUE
C
      IF (NPRT(9)) PRINT 1000
 1000 FORMAT(' *GNSLWD* ELASTIC SCATTERING')
C
      DO 24 J=4,9
      PV(J,1)=PV(J,MXGKPV)
 24   CONTINUE
      PV(10,1)=0.0
C
C --- VERY SIMPLE SIMULATION OF SCATTERING ANGLE AND ENERGY ---
C --- NONRELATIVISTIC APPROXIMATION WITH ISOTROPIC ANGULAR ---
C --- DISTRIBUTION IN THE CMS SYSTEM ---
  25  CALL GRNDM(RNDM,2)
      RAN=RNDM(1)
      COST1=-1.0+2.0*RAN
      EKA=1.0+2.0*COST1*ATNO2+ATNO2**2
      IF(EKA.LE.0.) GOTO 25
      COST=(ATNO2*COST1+1.0)/SQRT(EKA)
      IF (COST .LT. -1.0) COST=-1.0
      IF (COST .GT. 1.0) COST=1.0
      EKA=EKA/(1.0+ATNO2)**2
      EK=EK*EKA
      EN=EK+ABS(AMAS)
      P=SQRT(ABS(EN*EN-AMAS*AMAS))
      SINT=SQRT(ABS(1.0-COST*COST))
      PHI=RNDM(2)*TWPI
      PV(1,2)=SINT*SIN(PHI)
      PV(2,2)=SINT*COS(PHI)
      PV(3,2)=COST
      CALL DEFS1(2,MXGKPV,2)
      PU=SQRT(PV(1,2)**2+PV(2,2)**2+PV(3,2)**2)
      PX=PV(1,2)/PU
      PY=PV(2,2)/PU
      PZ=PV(3,2)/PU
      PV(1,1)=PX*P
      PV(2,1)=PY*P
      PV(3,1)=PZ*P
      PV(4,1)=EN
C
C --- STORE BACKSCATTERED PARTICLE FOR ATNO < 4.5 ---
      IF (ATNO2 .GT. 4.5) GO TO 27
C
      IF (NPRT(9)) PRINT 1001,ATNO2
 1001 FORMAT(' *GNSLWD* BACKSCATTERED PARTICLE STORED FOR ATNO ',G12.5)
C
      PV(1,2)=PV(1,MXGKPV)-PV(1,1)
      PV(2,2)=PV(2,MXGKPV)-PV(2,1)
      PV(3,2)=PV(3,MXGKPV)-PV(3,1)
      CALL LENGTX(2,PP)
      PV(9,2)=0.0
      PV(10,2)=0.0
      PV(7,2)=TOF
C
      IF (ATNO2 .GT. 3.5) GO TO 274
      IF (ATNO2 .GT. 2.5) GO TO 273
      IF (ATNO2 .GT. 1.5) GO TO 272
C
 271  CONTINUE
      PV(5,2)=RMASS(14)
      PV(4,2)=SQRT(PP*PP+PV(5,2)*PV(5,2))
      PV(6,2)=RCHARG(14)
      PV(8,2)=14.0
      GO TO 275
C
 272  CONTINUE
      PV(5,2)=RMASS(30)
      PV(4,2)=SQRT(PP*PP+PV(5,2)*PV(5,2))
      PV(6,2)=RCHARG(30)
      PV(8,2)=30.0
      GO TO 275
C
 273  CONTINUE
      PV(5,2)=RMASS(31)
      PV(4,2)=SQRT(PP*PP+PV(5,2)*PV(5,2))
      PV(6,2)=RCHARG(31)
      PV(8,2)=31.0
      GO TO 275
C
 274  CONTINUE
      PV(5,2)=RMASS(32)
      PV(4,2)=SQRT(PP*PP+PV(5,2)*PV(5,2))
      PV(6,2)=RCHARG(32)
      PV(8,2)=32.0
C
 275  CONTINUE
      INTCT=INTCT+1.0
      CALL SETCUR(1)
      NTK=NTK+1
      CALL SETTRK(2)
      GO TO 9999
C
C --- PUT QUANTITIES IN COMMON /RESULT/ ---
 27   CONTINUE
      IF (PV(10,1) .NE. 0.0) USERW=PV(10,1)
      SINL=PZ
      COSL=SQRT(ABS(1.0-SINL*SINL))
      IF (ABS(COSL) .LT. 1.E-10) GO TO 28
C
      SINP=PY/COSL
      COSP=PX/COSL
      GO TO 9999
C
 28   CONTINUE
      CALL GRNDM(RNDM,1)
      PHI=RNDM(1)*TWPI
      SINP=SIN(PHI)
      COSP=COS(PHI)
      GO TO 9999
C
C *** INELASTIC SCATTERING ON HEAVY NUCLEI ***
 29   CONTINUE
C
      IF (NPRT(9)) PRINT 1002
 1002 FORMAT(' *GNSLWD* INELASTIC SCATTERING ON HEAVY NUCLEUS')
C
C --- DECIDE BETWEEN SPALLATION OR SIMPLE NUCLEAR REACTION ---
      CALL GRNDM(RNDM,1)
      TEST1=RNDM(1)
      TEST2=4.5*(EK-0.01)
      IF (TEST1 .GT. TEST2) GO TO 40
C
C *** SPALLATION ***
C
      IF (NPRT(9)) PRINT 1003
 1003 FORMAT(' *GNSLWD* SPALLATION')
C
      PV( 1,MXGKPV)=P*PX
      PV( 2,MXGKPV)=P*PY
      PV( 3,MXGKPV)=P*PZ
      PV( 4,MXGKPV)=EN
      PV( 5,MXGKPV)=AMAS
      PV( 6,MXGKPV)=NCH
      PV( 7,MXGKPV)=TOF
      PV( 8,MXGKPV)=IPART
      PV( 9,MXGKPV)=0.0
      PV(10,MXGKPV)=USERW
C
C --- FERMI-MOTION AND EVAPORATION ---
      TKIN=CINEMA(EK)
      ENP(5)=EK+TKIN
C --- CHECK FOR LOWERBOUND OF EKIN IN CROSS-SECTION TABLES ---
      IF (ENP(5) .LE. TEKLOW) ENP(5)=TEKLOW
      ENP(6)=ENP(5)+ABS(AMAS)
      ENP(7)=ENP(6)*ENP(6)-AMASQ
      ENP(7)=SQRT(ENP(7))
      TKIN=FERMI(ENP(5))
      ENP(5)=ENP(5)+TKIN
C --- CHECK FOR LOWERBOUND OF EKIN IN CROSS-SECTION TABLES ---
      IF (ENP(5) .LE. TEKLOW) ENP(5)=TEKLOW
      ENP(6)=ENP(5)+ABS(AMAS)
      ENP(7)=ENP(6)*ENP(6)-AMASQ
      ENP(7)=SQRT(ENP(7))
      TKIN=EXNU(ENP(5))
      ENP(5)=ENP(5)-TKIN
C --- CHECK FOR LOWERBOUND OF EKIN IN CROSS-SECTION TABLES ---
      IF (ENP(5) .LE. TEKLOW) ENP(5)=TEKLOW
      ENP(6)=ENP(5)+ABS(AMAS)
      ENP(7)=ENP(6)*ENP(6)-AMASQ
      ENP(7)=SQRT(ENP(7))
C
C --- NEUTRON CASCADE ---
      K=2
      CALL VZERO(IPA(1),MXGKCU)
      CALL CASN(K,INT,NFL)
      GO TO 9999
C
 40   CONTINUE
      IF (NPRT(9)) PRINT 1004
 1004 FORMAT(' *GNSLWD* NUCLEAR REACTION')
      CALL NUCREC(NOPT,1)
      IF (NOPT .NE. 0) NUCFLG=1
      IF (NOPT .EQ. 0) GO TO 30
C
 9999 CONTINUE
      END
